---
title: "Development and characterization of microsatellite markers for population
  Gracillaridae)"
author: "Marynold Purificacion1#, Roslina Binti Mohd Shah2#, Thierry De Meeûs3, Saripah
  Binti Bakar4, Anisah Bintil Savantil2, Meriam Mohd Yusof2, Divina Amalin1, Hien
  Nguyen5, Endang Sulistyowati6, Aris Budiman6, Arni Ekayanti7, Jerome Niogret8, Sophie
  Ravel3, Marc J.B. Vreysen9, Adly M.M. Abd-Alla9"
date: "2023-02-23"
output: word_document
---

```{}
setwd("C:/Users/abdallaa/OneDrive - IAEA/My_passport_6/Marynold/analysis")

library(ggplot2)
library(MASS)
library(rmarkdown)
library(knitr)
library(lme4)
library(MuMIn)
library(ggthemes) # Load
library(datasets)
library(plyr)
library(dplyr)
library(tidyverse)
library(geosphere)

```
## Analysis CPB 11 loci; ajusted p value for LD

```{r }
p.adjust(c(0.0125, 0.5208, 0.4659, 0.1141, 0.2034, 0.6403, 0.4994, 0.633,	0.8052,	0.4462,	0.0002,	0.8846,	0.5612, 0.0445,	0.2988,	0.9458,	0.1942, 0.0614,	0.2284, 0.7694,	0.4195,	0.3797,	0.3607,	0.4739,	0.4429,	0.2802,	0.2529,	0.1332,	0.4537,	0.1205,	0.3059,	0.0323,	0.6229,	0.337, 1, 0.1831, 0.5724, 1, 0.159, 0.7901, 0.5994, 0.3967, 0.3837, 0.1554, 0.0971, 0.1102, 0.0876, 0.2744, 0.3586, 0.0002, 0.8103, 0.0031, 0.0001, 0.4601, 0.3488), method="BY")

```
## Analysis Wahlund effect test for CPB with 11 loci

```{r }

# Wahlund effect test for CPB with 11 loci
wahlund <- read.csv("CPB_11loci_Wahlund.csv",sep=",", row.names=NULL)
wahlund
with(wahlund, cor.test(Ht, nLDsig, alternative="less", method="spearman"))

```

## Analysis CPB 11 loci; ajusted p value for observed heterozygocity

```{r }
p.adjust(c(0.3627, 0.9999, 0.924, 1, 0.9971, 0.9987, 0.9935, 1, 0.1311, 0.9926, 0.9934), method="BH")

```

## Analysis the correlation betweel FIS and FST and number of blanks for CPB 11 loci and Fis and number of blanks with 7 loci

```{r }

# Correlation between Fis and fst for cPB with 11 loci

fisfst <- read.csv("cpb_11loci_fis_fst.csv",sep=",", row.names=NULL)
fisfst
with(fisfst, cor.test(FIS, FST, alternative="greater", method="spearman"))


# Correlation between Fis and number of blank for cPB with 11 loci

fisblank11 <- read.csv("Fis_blanks_11loci.csv",sep=",", row.names=NULL)
fisblank11
with(fisblank11, cor.test(FIS, Blanks, alternative="greater", method="spearman"))

# Correlation between Fis and number of blank for cPB with 7 loci
fisblank7 <- read.csv("Fis_blanks_7loci.csv",sep=",", row.names=NULL)
fisblank7
with(fisblank7, cor.test(FIS, Blanks, alternative="greater", method="spearman"))

```

## Analysis for short allel dominance (SAD) for CPB with 11 loci


```{r }
SAD <- read.csv("SAD.csv",sep=",", row.names=NULL)
SAD
SAD$weight<- as.numeric(SAD$weight)
SAD$Allele <- as.numeric(SAD$Allele)

SAD_Cpb14 <- subset(SAD,loci=="Cpb14")
SAD_Cpb14
with(SAD_Cpb14, cor.test(Allele, Capf, alternative="less", method="spearman"))


SAD_Cpb55 <- subset(SAD,loci=="Cpb55")
SAD_Cpb55
with(SAD_Cpb55, cor.test(Allele, Capf, alternative="less", method="spearman"))

SAD_Cpb54 <- subset(SAD,loci=="Cpb54")
SAD_Cpb54
with(SAD_Cpb54, cor.test(Allele, Capf, alternative="less", method="spearman"))
lm_SAD_54 <- lm(Smallf ~ Allele,	 data=SAD_Cpb54, weights=weight)
summary(lm_SAD_54)

SAD_Cpb52 <- subset(SAD,loci=="Cpb52")
SAD_Cpb52
with(SAD_Cpb52, cor.test(Allele, Capf, alternative="less", method="spearman"))

SAD_Cpb84 <- subset(SAD,loci=="Cpb84")
SAD_Cpb84
with(SAD_Cpb84, cor.test(Allele, Capf, alternative="less", method="spearman"))

SAD_Cpb84 <- subset(SAD,loci=="Cpb84")
SAD_Cpb84
with(SAD_Cpb84, cor.test(Allele, Capf, alternative="less", method="spearman"))
lm_SAD_84 <- lm(Smallf ~ Allele,	 data=SAD_Cpb84, weights=weight)
summary(lm_SAD_84)


SAD_Cpb122 <- subset(SAD,loci=="Cpb122")
SAD_Cpb122
with(SAD_Cpb122, cor.test(Allele, Capf, alternative="less", method="spearman"))

SAD_Cpb135 <- subset(SAD,loci=="Cpb135")
SAD_Cpb135
with(SAD_Cpb135, cor.test(Allele, Capf, alternative="less", method="spearman"))
lm_SAD_135 <- lm(Smallf ~ Allele,	 data=SAD_Cpb135, weights=weight)
summary(lm_SAD_135)

SAD_Cpb160 <- subset(SAD,loci=="Cpb160")
SAD_Cpb160
with(SAD_Cpb160, cor.test(Allele, Capf, alternative="less", method="spearman"))
lm_SAD_160 <- lm(Smallf ~ Allele,	 data=SAD_Cpb160, weights=weight)
summary(lm_SAD_160)


SAD_Cpb139 <- subset(SAD,loci=="Cpb139")
SAD_Cpb139


```

## Analysis the bionomial test of the expected blank for null alleles of CPB with 11 loci 

```{r }
stut <- read.csv("cpb_11loci_stutt_bionom.csv",sep=",", row.names=NULL)
stut
stut_cpb14<-subset(stut,Loci=="Cpb14")
stut_cpb14
with(stut_cpb14,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb55<-subset(stut,Loci=="Cpb55")
stut_cpb55
with(stut_cpb55,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb54<-subset(stut,Loci=="Cpb54")
stut_cpb54
with(stut_cpb54,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb62<-subset(stut,Loci=="Cpb62")
stut_cpb62
with(stut_cpb62,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb84<-subset(stut,Loci=="Cpb84")
stut_cpb84
with(stut_cpb84,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb122<-subset(stut,Loci=="Cpb122")
stut_cpb122
with(stut_cpb122,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb112<-subset(stut,Loci=="Cpb112")
stut_cpb112
with(stut_cpb112,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb135<-subset(stut,Loci=="Cpb135")
stut_cpb135
with(stut_cpb135,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb160<-subset(stut,Loci=="Cpb160")
stut_cpb160
with(stut_cpb160,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb139<-subset(stut,Loci=="Cpb139")
stut_cpb139
with(stut_cpb139,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

stut_cpb190<-subset(stut,Loci=="Cpb190")
stut_cpb190
with(stut_cpb190,binom.test(ObsBlanks, N, f.ExpBlanks., alternative="less"))

# adjusted P value
p.adjust(c(0.3566,0.7657,0.9981,0.9985,0.9454,0.9968,1,0.02423,0.3383,0.9968,0.5842), method="BH")

```

## Analysis the correlation between GST and Hs for CPB with 7 loci

```{r }

#analysis the correlation between GST and HS for CPB with 7 loci

GstHs7 <- read.csv("Gst_Hs_7loci.csv",sep=",", row.names=NULL)
GstHs7
with(GstHs7, cor.test(GST, HS, alternative="less", method="spearman"))

```
## Analysis the correlation between FIS and FST  and number of blanks for CPB with 7 loci

```{r }

#correlation between FIS and FST for CPB with 7 loci

fisfst7 <- read.csv("cpb_7loci_fis_fst.csv",sep=",", row.names=NULL)
fisfst7
with(fisfst7, cor.test(FIS, FST, alternative="greater", method="spearman"))

#correlation between FIS and number of blanks with 5 loci
fisblank5 <- read.csv("Fis_blanks_5loci.csv",sep=",", row.names=NULL)
fisblank5
with(fisblank5, cor.test(FIS, Blanks, alternative="greater", method="spearman"))

```
## Analysis the Geographical discatnce between location for CPB with 7 loci

```{r }
library(geosphere)
Longlat11<-read.table("Longlat11.txt", header = TRUE)
Longlat22<-read.table("Longlat22.txt", header = TRUE)
tabdistgeo<-data.frame(distGeo(Longlat11, Longlat22))
write.table(tabdistgeo,"tabdistgeo.txt",col=NA,sep="\t",dec=".")
dgeo <- read.table("tabdistgeo.txt",header = TRUE)
dgeo

```

## Analysis  adjusted p value for LFB for LD values and for observed heterogeneity

```{r }
p.adjust(c(00.8517, 0.5328, 0.3333, 0.6438, 0.7802, 0.7005, 0.9416, 0.5651, 0.4517, 0.8289, 0.0069, 1, 0.8056, 0.691, 0.6992, 1, 0.0503, 0.3089, 0.6856, 1, 0.6093, 0.7343, 1, 0.7965, 0.1931, 1, 0.8441, 0.7123, 1, 0.81, 1, 1, 0.574, 0.5147, 1, 1, 0.1716, 1, 0.4394, 0.5623, 1, 0.5261, 1, 0.5381, 0.6493, 0.8413, 1, 1, 1, 0.5727, 0.9562, 0.4851, 0.7081, 0.8287, 1), method="BY")

# adjusted P value for observed heterogeneity
p.adjust(c(0.7664, 0.2309, 0.2371, 1, 0.1652, 0.6948, 0.1608, 1, 0.01838, 0.9332, 0.07359), method="BH")

```


## Analysis for short allel dominance (SAD) for LFB with 11 loci

```{r }
LFB_SAD <- read.csv("LFB_SAD.csv",sep=",", row.names=NULL)
LFB_SAD
LFB_SAD$weight<- as.numeric(LFB_SAD$weight)
LFB_SAD$Allele <- as.numeric(LFB_SAD$Allele)

LFB_SAD_Cpb14 <- subset(LFB_SAD,loci=="Cpb14")
LFB_SAD_Cpb14
lm_LFB_SAD_Cpb14 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb14, weights=weight)
summary(lm_LFB_SAD_Cpb14)


LFB_SAD_Cpb55 <- subset(LFB_SAD,loci=="Cpb55")
LFB_SAD_Cpb55
lm_LFB_SAD_Cpb55 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb55, weights=weight)
summary(lm_LFB_SAD_Cpb55)

LFB_SAD_Cpb54 <- subset(LFB_SAD,loci=="Cpb54")
LFB_SAD_Cpb54
lm_LFB_SAD_Cpb54 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb54, weights=weight)
summary(lm_LFB_SAD_Cpb54)

LFB_SAD_Cpb62 <- subset(LFB_SAD,loci=="Cpb62")
LFB_SAD_Cpb62
lm_LFB_SAD_Cpb62 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb62, weights=weight)
summary(lm_LFB_SAD_Cpb62)

LFB_SAD_Cpb84 <- subset(LFB_SAD,loci=="Cpb84")
LFB_SAD_Cpb84
lm_LFB_SAD_Cpb84 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb84, weights=weight)
summary(lm_LFB_SAD_Cpb84)

LFB_SAD_Cpb122 <- subset(LFB_SAD,loci=="Cpb122")
LFB_SAD_Cpb122
lm_LFB_SAD_Cpb122 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb122, weights=weight)
summary(lm_LFB_SAD_Cpb122)


LFB_SAD_Cpb112 <- subset(LFB_SAD,loci=="Cpb112")
LFB_SAD_Cpb112
lm_LFB_SAD_Cpb112 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb112, weights=weight)
summary(lm_LFB_SAD_Cpb112)

LFB_SAD_Cpb135 <- subset(LFB_SAD,loci=="Cpb135")
LFB_SAD_Cpb135
lm_LFB_SAD_Cpb135 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb135, weights=weight)
summary(lm_LFB_SAD_Cpb135)


LFB_SAD_Cpb160 <- subset(LFB_SAD,loci=="Cpb160")
LFB_SAD_Cpb160
lm_LFB_SAD_Cpb160 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb160, weights=weight)
summary(lm_LFB_SAD_Cpb160)


LFB_SAD_Cpb139 <- subset(LFB_SAD,loci=="Cpb139")
LFB_SAD_Cpb139
lm_LFB_SAD_Cpb139 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb139, weights=weight)
summary(lm_LFB_SAD_Cpb139)

LFB_SAD_Cpb190 <- subset(LFB_SAD,loci=="Cpb190")
LFB_SAD_Cpb190
lm_LFB_SAD_Cpb190 <- lm(Smallf ~ Allele, data=LFB_SAD_Cpb190, weights=weight)
summary(lm_LFB_SAD_Cpb190)

```


## Analysis the correlation of FIS and number of blanks for LFB with 11 loci



```{r }
# correlation of FIS and number of blanks for LFB with 11 loci

lfbfisblank11 <- read.csv("LFB_11loci_fis_blank.csv",sep=",", row.names=NULL)
lfbfisblank11
with(lfbfisblank11, cor.test(FIS, Blanks, alternative="greater", method="spearman"))


```


## Analysis of male female divergence for supplementary file 2


```{r }
#male female divergence_supplementary file 2
subfile2 <- read.csv("malefemale_divergence.csv",sep=",", row.names=NULL)
subfile2
subfile2_Cpb14 <- subset(subfile2,loci=="Cpb14")
subfile2_Cpb14
with(subfile2_Cpb14, wilcox.test(FIS.Females, FIS.Males, alternative='two.sided', paired=TRUE))

subfile2_Cpb112 <- subset(subfile2,loci=="Cpb112")
subfile2_Cpb112
with(subfile2_Cpb112, wilcox.test(FIS.Females, FIS.Males, alternative='two.sided', paired=TRUE))

# Bh adjusted P value for supplementary file 2
p.adjust(c(0.6874, 0.6874, 0.3165, 0.6874, 0.2241, 0.0151, 0.0684, 0.6874, 0.179, 0.6874, 0.2808), method="BH")
subfile2_Cpb122 <- subset(subfile2,loci=="Cpb122")
subfile2_Cpb122
with(subfile2_Cpb122, wilcox.test(FIS.Females, FIS.Males, alternative='two.sided', paired=TRUE))


```


## 



```{r }


```

